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Electron spin relaxation due to the D'yakonov-Perel' mechanism is investigated in bilayer graphene 
with only the lowest conduction band being relevant. The spin-orbit coupling is constructed from 
the symmetry group analysis with the parameters obtained by fitting to the numerical calculation 
according to the latest report by Konschuh et al. [Phys. Rev. B 85, 115423 (2012)] from first 
principles. In contrast to single-layer graphene, the leading term of the out-of-plane component of the 
spin-orbit coupling in bilayer graphene shows a Zeeman-like term with opposite effective magnetic 
fields in the two valleys. This Zeeman-like term opens a spin relaxation channel in the presence of 
intervalley scattering. It is shown that the intervalley electron-phonon scattering, which has not 
been reported in the previous literature, strongly suppresses the in-plane spin relaxation time at high 
temperature whereas the intervalley short-range scattering plays an important role in the in-plane 
spin relaxation especially at low temperature. A marked nonmonotonic dependence of the in-plane 
spin relaxation time on temperature with a minimum of several hundred picoseconds is predicted in 
the absence of the short-range scatterers. This minimum is comparable to the experimental data. 
Moreover, a peak in the electron density dependence of the in-plane spin relaxation time at low 
temperature, which is very different from the one in semiconductors, is predicted. We also find a 
rapid decrease in the in-plane spin relaxation time with increasing initial spin polarization at low 
temperature, which is opposite to the situation in both semiconductors and single-layer graphene. 
A strong anisotropy between the out-of- and in-plane spin relaxations at high temperature is also 
revealed with the out-of-plane spin relaxation time being about two orders of magnitude larger than 
the in-plane one. Detailed comparisons of the temperature and electron density dependences of the 
spin relaxation with the existing experiments [Phys. Rev. Lett. 107, 047206 (2011), Phys. Rev. 
Lett. 107, 047207 (2011) and Nano Lett. 11, 2363 (2011)] are reported. 

PACS numbers: 72.25.Rb, 81.05.ue, 71.10.-w, 71.70.Ej 



I. INTRODUCTION 

In recent years, single-layer graphene has received 
much attention due to its two dimensionality, Dirac- 
like energy spectrum and outstanding spin-coherence 
properties.— ~— Specifically, the hyperfine interaction and 
spin-orbit coupling (SOC) are considerably weakj 1 ^— 
suggesting a long spin relaxation time (SRT) in single- 
layer graphene. This long SRT makes single-layer 
graphene a promising candidate for spintronic applica- 
tion. So far, a number of experiments on spin relax- 
ation have been carried out, with the SRTs being orders 
of magnitude shorter than theoretical predictions^ - — 
This pronounced discrepancy has provoked many in- 
vestigations on the extrinsic effects of spin relaxation 
such as adatomS ) 24 ' 42 " — curvature ) 16 i 48 ~ — substrate 
effects^iiS 2 . and contacts 1 4 ^ Consisting of two 
layer graphene, bilayer graphene (BLG) also possesses 
good spin-coherence propertiesj 20 ' 21 i 54 ~ — Additionally, in 
contrast to single-layer graphene, BLG has a tunable gap 
caused by an out-of-plane electric field ) 55 ' 56 i 58 ~ — which 
can be used to turn on or off the electronic transport. In 
addition, in the absence of space inversion symmetry, the 
intrinsic SOC in the pseudospin space of BLG provides an 
out-of-plane effective magnetic field near the Dirac points 
whereas the intrinsic SOC of single-layer graphene only 
induces a shift of the energy spectrum i 36 ' 44 i 55 Moreover, 



the electron-phonon coupling becomes more complex in 
BLG due to the existence of more phonon modest - — 
So far, some efforts have been made to investigate the in- 
travalley electron-phonon coupling i^ 3 - - — Specifically, the 
intravalley electron-acoustic (AC) phonon coupling has 
been calculated by Viljas and Heikkila 63 with the con- 
tinuum theory and also by Borysenko et al£± from first 
principles. To explore the intravalley electron-optical 
(OP) phonon coupling, the tight-binding model, density 
functional theory and first principles approach are em- 
ployed by Viljas and Heikkila,,- 3 Cappelluti et al£^ and 
Borysenko et al.^ respectively. However, the intervalley 
electron-phonon coupling, which has been shown to play 
an important role in spin relaxation in rippled single- 
layer graphene^ has not been reported in the previous 
literature. All these intriguing features make BLG an 
attractive material for investigation. 

Very recently, some attention has been devoted to the 
spin relaxation in BLGi 3 ^ 39 ' 66 " — Experimentally, the 
SRTs of BLG on Si0 2 substrate^L 3 ^ or j n f ree ly sus- 
pended BLG^ 8 " are reported to be of the order of 0.01- 
1 ns. Both Yang et al£^ and Han and Kawakami 3 - 7 - mea- 
sured the temperature and carrier dependences of the 
SRT. Yang et al. 66 showed that the SRT is weakly de- 
pendent on temperature. However, Han and Kawakami 37 
reported a marked decrease of the SRT with the increase 
of the temperature. As for the carrier density depen- 
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dence, the decrease of the SRT with increasing carrier 
density at low temperature is reported in both exper- 
iments. However, at room temperature, Yang et a/.— 
observed an increase of the SRT with increasing density 
whereas the SRT measured by Han and Kawakami^I ex- 
hibits a marginal density dependence. In addition, Yang 
et al. 6e also showed that the SRT scales inversely with the 
mobility at both room and low temperatures. In contrast 
to the above two experiments based on the micromechan- 
ically exfoliated BLG samples, Avsar et a/.— carried out 
the experiment in chemical vapor deposition grown BLG. 
The carrier density and temperature dependences of the 
SRT measured are comparable to the ones by Yang et 
al£& both quantitatively and qualitatively. Motivated 
by these experiments, Diez and Burkard^ calculated 
the electron spin relaxation due to the D'yakonov-Perel' 
(DP) mechanism 69 with only the long-range electron- 
impurity scattering included in BLG. They claimed to 
obtain a good agreement with the experimental data of 
Yang et al.£& It is noted that the effective SOC (re- 
ferred to as the spin-orbit field in Ref. l67h used in their 
calculation is induced by an interlayer bias voltage to- 
gether with the intrinsic SOC in the pseudospin spaced 4 
However, the intrinsic SOC in the pseudospin space be- 
comes incomplete in the presence of the interlayer bias 
since the interlayer bias breaks the space inversion sym- 
metry and leads to additional extrinsic SOC terms in 
the pseudospin space. 55 Also even with only the intrinsic 
SOC in the pseudospin space in their calculation, they 
overlooked the out-of-plane component of the effective 
SOC. As will be shown later, the leading term of this 
component leads to a Zeeman-like term in the two val- 
leys and makes a dominant contribution to the in-plane 
spin relaxation. Moreover, they neglected the electron- 
phonon and electron-electron Coulomb scatterings. In 
fact, the electron-electron Coulomb scattering has been 
demonstrated to be very important for the spin relax- 
ation in semiconductor systems^ - - — and also single-layer 
graphene at high temperature^ 6 - The electron-phonon 
scattering has also been shown to play an important role 
in spin relaxation in both semiconductors^ - and single- 
layer graphene^ 6 - at high temperature. Therefore, both 
have to be included in studying the spin relaxation. 



In the present work, with the electron-electron 
Coulomb, (both the intra- and inter-valley) electron- 
phonon, possible short-range 9 - as well as long-range 
electron-impurity scatterings explicitly included, we in- 
vestigate the electron spin relaxation due to the DP 
mechanism in BLG with only the lowest conduction band 
being relevant by the kinetic spin Bloch equation (KSBE) 
approach. 70 The SOC of the lowest conduction band near 
the Dirac points f2 M (k) is constructed from the symme- 
try group C3 in the presence of an out-of-plane electric 
fields Its three components f2f(k) (i — x,y,z) can be 



written as (see Appendix |A")) 

fl£(k) = ai(/c)sin6> k + ^[a 2 (/c)sin2(9 k 

+ a 3 (*)sin40 k ], (1) 

fi£(k) = — ai(fc) cos6> k + ^[aa(fc) cos2# k 

-a 3 (fc)cos40 k ], (2) 

n{f(k) = Mk)+ /3 2 (£Ocos30 k , (3) 

with jj, = 1(— 1) standing for K(K') valleys and k and 
k representing the magnitude and the polar angle of 
the momentum k relative to K(K') points, respectively. 
The coefficients a,(fc) (i = 1-3) and /3i.2 (^) are given 
in Appendix [A] by fitting to the numerical results calcu- 
lated according to the latest report by Konschuh et a/.— 
from first principles with both the intrinsic and extrin- 
sic SOC terms in the pseudospin space included. The 
out-of-plane component n^(k) [Eq. which has not 
been reported in the previous literature, is induced by 
the intrinsic SOC in the pseudospin space in the absence 
of space inversion symmetry, whereas the in-plane com- 
ponent is contributed by both the intrinsic and extrinsic 
SOCs in the pseudospin space. This is very different 
from the case in single-layer graphene where the intrin- 
sic SOC in the pseudospin space only induces a shift 
of the energy spectrum in the conduction band i 36 i 44 i 55 
We find that the magnitudes of both the out-of- and in- 
plane components decrease with the increase of k at large 
momentum, indicating a suppression of the inhomoge- 
neous broadening 71 by increasing momentum. 55 This is 
very different from the case in both semiconductors and 
single-layer graphene i 36 i 70 Additionally, the leading term 
of the out-of-plane component, i.e., /i/3i(fc), serves as a 
Zeeman-like term with opposite effective magnetic fields 
perpendicular to BLG plane in the two valleys, similar to 
the case in rippled single-layer graphene^ - This Zeeman- 
like term, together with the intervalley scattering, opens 
a spin relaxation channel. This spin relaxation channel 
suppresses the in-plane SRT significantly at high temper- 
ature by the intervalley electron-phonon scattering and 
also at low temperature by the short-range scattering. 9 
In addition, we find that although such Zeeman-like 
terms also exist in the in-plane components of the SOC, 
their contributions to the out-of-plane spin relaxation are 
marginal since these terms are of the high order of the 
momentum. This leads to a strong anisotropy between 
the in-plane and out-of-plane spin relaxations. We also 
find that the in-plane SRT shows a marked nonmonotonic 
dependence on temperature with a minimum down to 
several hundred picoseconds in the absence of short-range 
scatterers. This minimum is comparable to the experi- 
mental data. The nonmonotonic behavior results from 
the crossover between the weak and strong intervalley 
electron-phonon scattering. Moreover, we predict a peak 
in the density dependence at low temperature, which is 
very different from the one in semiconductors 1 77 i 79 As for 
high temperature, the SRT shows a monotonic increase 
with increasing density. We also find that the in-plane 
SRT decreases rapidly with increasing initial spin polar- 
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ization at low temperature, which is very different from 
the previous studies in both semiconductor s 72 ' 76 ' 80 ' 81 and 
single-layer graphene, 36 Finally, we report detailed com- 
parisons of the temperature and electron density depen- 
dences with the existing experiments ! 37 ' 38 ' 66 

This paper is organized as follows. In Sec. II, we 
present our model and construct the KSBEs. The main 
results are given in Sec. III. In Sec. IIIA, we investigate 
the temperature, electron density and also initial spin 
polarization dependences of the electron spin relaxation 
in the absence of short-range scatterers. The anisotropy 
of the spin relaxation is also addressed in this part. In 
Sec. IIIB, we show detailed comparisons of the tempera- 
ture and electron density dependences of the spin relax- 
ation with the existing experiments with the inclusion of 
short-range scatterers. We summarize in Sec. IV. 



II. MODEL AND KSBES 

We start our investigation from the AB-stacked BLG 
on SiC>2 substrate. The spinless 7r-band structure can be 
described by the effective 4x4 Hamiltonian: 55 
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in the on-site orbital Bloch basis (k), "Jbi (k), *S?a 2 (k) 
and 1 $b 2 (k), where Ai refers to the A sublattice in the 
lower layer, B2 to the B sublattice in the upper layer, 
etc., and k represents the two-dimensional wave vec- 
tor counted from the T point. Here, 70 and 71 de- 
scribe the nearest-neighbor intralayer and interlayer hop- 
pings whereas 73 and 74 stand for the indirect hop- 
pings between two layers. A represents the asymme- 
try in the energy shift of the on-site energy due to 
the interlayer hopping. The potential V — eE z d e g/2 
with e, E z and d c g representing the electron charge 
(e > 0), the out-of-plane electric field and the effec- 
tive electrostatic bilayer distance, respectively, /(k) = 
e l ik ky \i 4- 2e~ J ^f i * ; y cos(#fc x )] is the nearest-neighbor 
structural function of the graphene hexagonal lattice with 
a denoting the lattice constant. Near the Dirac points, 
to the first order of the momentum k originated from the 
K(K') points, /(k) = p(k) = -V3a(pk x - ik v )/2. Then 
the effective 4x4 Hamiltonian near the Dirac points can 
be written as 
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A + V 70P 74P* 71 

loP* +V 73P 74P* 

J4P 13P* -V 70P 

7i 1+P loP* A-V 



(5) 



By exactly diagonalizing this Hamiltonian, one obtains 
the eigenvalue e M „k and the eigenfunction ip£ v with v be- 
ing the band index including two conduction bands and 



two valence ones. In our calculation, we concentrate on 
the heavily n-doped case with only the lowest conduction 
band being relevant whereas other bands are separated 
by large energy intervals. For convenience, the eigen- 
value and the corresponding eigenfunction of the lowest 
conduction band at k are denoted by e M k and t/£ , respec- 
tively. 

By incorporating the spin degree of freedom, the effec- 
tive Hamiltonian of the lowest conduction band is given 
by 



cn 



e Alk + n^(k)- ( r/2. 



(6) 



Here, <r are the Pauli matrices. i~2 M (k) stands for 
the effective magnetic field of the SOC, which is con- 
structed from the symmetry analysis with the parame- 
ters obtained approximately by fitting to numerical re- 
sults based on the report by Konschuh et al.£^ The an- 
alytical form of the SOC with three components fif (k) 
(i = x, y, z) is shown in Eqs. ([T][3]) and the explicit numer- 
ical calculation is given in Appendix [AJ It is noted that 
at large momentum, the analytical result agrees with the 
numerical one fairly well (see Appendix [AT) . 

We then construct the microscopic KSBEsZ - to study 
the electron spin relaxation in BLG. The KSBEs read 70 



(7) 



where p^ represent the density matrices of electrons 
with the diagonal terms p^.aa = f^a (c = ±5) de- 
scribing the distribution functions and the off-diagonal 
ones /Vk,(i/2)(-i/2) = p£k,(-i/2)(i/2) denoting the spin 
coherence. The coherent term is given by 



d t p^ 



coh 



(8) 



where [A, B] = AB — BA is the commutator; 
£>k = _ Sk' ^klk'^kk'^k' represents the Coulomb 
Hartree-Fock (HF) termi 7 - 2 - V£ k , denotes the screened 
Coulomb potential with its form given in Appendix IB"1 

The form factor i£ k , = IV'kVk-l 2 - 9tP/<k|scat 
are the scattering terms including the electron- 
electron Coulomb (IT^V/ | 2 ), long-range electron- impurity 
(l^kU'l 2 )' intravalley electron-AC phonon (|M^ , k , | 2 ), 
electron-OP phonon (|M"° k P M , k , | 2 ), electron-remote- 
interfacial (RI) phonon (|M^ k , k , | 2 ) and especially 
the intervalley electron-phonon scattering including 

the electron-KA^ phonon (|M^, k ,| 2 ), electron-KA 2 

phonon (|M*^, k ,| 2 ), electron-KE' phonon (|M™;, k ,| 2 ) 
and electron-KZ phonon (\M^ , k , | 2 ) scatterings. Their 
detailed expressions can be found in Ref. I36U 2 - Here, 
KA[ , KA 2 and KE' are the in-plane phonon modes corre- 
sponding to the representations A[, A' 2 and E' of group 
Z?3/j, respectively whereas KZ stands for the out-of-plane 
phonon mode i 64 ' 83 " — It is noted that the intravalley 
electron-phonon scattering has been reported in the pre- 
vious literature as mentioned in the introduction ! 63 ' 
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However, to the best of our knowledge, there is no spe- 
cific report to the intervalley electron-phonon scattering. 
In this work, we derive the intervalley electron-phonon 
scattering matrix elements using the tight-binding model 
according to the arXiv version of Ref. 163J. The scattering 
matrix elements are given by 

m 2 „ 



M 



KA[ , 
//k^'k' I 



Apfl 
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|< t (2v / 37o4 3 +«7;7SA 4 ) 



\M 
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(11) 



7i 



P 2i1kZ 



(12) 

where p is the mass density of the BLGjS 3 " sIka' , ^ka' j 
SIke' and SIkz are the energy spectra of KA\, KAj, 
KE', and KZ phonon modes, respectively) 64 ' 83 ' 84 7,' (i = 
0, 1,3,4) is the derivative of 7; with respect to the cor- 
responding hopping bond length^ I4 is the bond length 
of the inter layer hopping 74; 7^ (i — 0-3, 5), Oq 1 and 
tTp 3 are 4x4 Dirac matrices given in Appendix [Cj 86 We 
will show that the above intervalley electron-phonon scat- 
tering plays a significant role in the in-plane spin relax- 
ation in our investigation. The remaining scattering ma- 



trix elements such as \V^ k , | 2 , 



|[/£ k ,| 2 , and |M^, k , 



(A = AC, OP and PJ) are given in Appendix [B] It is 
noted that we also include the short-range scattering in 
our calculation with the scattering term laid out later in 
this paper. 



TABLE I: Parameters used in the calculation (from Ref. [55 
unless otherwise specified). 



A 


9.7 meV 


70 


2.6 eV 


n 


0.339 eV 


7 3 


0.28 eV 


74 


-0.14 eV 


An 


12 peV 


A12 


10 peV 


Ao 


(5 + 10.E*) peV 


a;, 


(5-10E*) peV 


A:, 


1.5E* peV 


A 4 


(-12-3SJ) /xeV A 4 


(-12 + 3BJ) peV 


a 


2.46 A 




0.1 nm 


d 


0.4 nm a 


r s 


0.8 a 




8 x 10 5 m/s b 


Vph 


2 x 10 4 m/s a 


D AC 


15 eV c 


P 


1.52 x 10~ 7 g/cm 




59 meV 




155 meV a 


Hi 


5.4 x 10" 3 a 


92 


3.5 x 10" 2 a 


ft LT 


196 meV a 




104 meV c 


7o 


4.4 eV/A e 


ix 


0.61 eV/A e 


73 


0.54 eV/A e 


f 

74 


0.3 eV/A e 




161.2 meV 




125 meV c 


s7kE' 


152 meV c 




65 meV c 


h 


3.64 A 


u 


3.64 A 


z, 


1 







a Reference 36 . 

b Reference 87 . 

c Reference 64 . 
Reference 63 . 

c Reference 65 . 
It is noted that E* 
unit of V/nm. 



is the magnitude of E z in the 



A. Spin relaxation in BLG with high mobilities 

In our calculation, the out-of-plane electric field is 
chosen to be a typical value E z — 0.14 V/nrr*2I and 
the long-range impurity density is taken to be Ni — 
1.5x 10 11 cm~ 2 . The corresponding mobility in our inves- 
tigation is of the order of 10 4 cm 2 /(Vs), which is about 
1-2 orders of magnitude larger than those reported in 
the existing experiments in BLG on Si02 substrat o 37 ' 38 ' 66 
and comparable to the ones in the experiments in freely 
suspended BLGj££ It is noted that the short-range scat- 
tering is not included in this subsection. 



III. NUMERICAL RESULTS 



By numerically solving the KSBEs, one obtains the 
time evolution of spin polarization along direction n 
p (t) = E M k Tr[^ k (t)<r • n]/N e with N e being the elec- 
tron density. Then the SRT can be determined from the 
slope of the envelope of the spin polarization P(t)JS. The 
initial spin polarization P(Q) is set to be 2.5 % and the 
spin-polarization direction n is along the x-axis unless 
otherwise specified. It is noted that all the parameters 
used in our calculation are listed in Table I. It is also 
noted that the SOC is calculated both fully numerically 
following the approach by Konschuh et al£^ and analyt- 
ically by using Eqs. ([US]) . The SRT calculated with the 
analytical form of the SOC and with the explicit numer- 
ical one agrees fairly well with each other. 



1. Temperature dependence of spin relaxation 

We first investigate the temperature dependence of the 
spin relaxation. In Fig. [1] the SRT along the x-axis r Sx 
is plotted as function of temperature T. It is seen that 
T,, m changes little with the temperature at low temper- 
ature (i.e., T < 60 K). However, as T becomes larger, 
the SRT decreases dramatically with increasing temper- 
ature. In particular, the SRT reaches down to several 
hundred picoseconds at T — 200 K, which is comparable 
to the experimental values i 37 ' 38 ' 66 ' 68 A detailed compar- 
ison with the existing experiments will be discussed in 
Sec. IIIB. When T further increases, the SRT begins to 
increase with T. 
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scattering time, respectively. According to the report by 
Zhang et al.^ the in-plane SRT is given by 



FIG. 1: (Color online) Total SRT along the x-axis t Sx (x) and 
the calculation without the long-range electron-impurity (•), 
electron-AC-phonon (a), electron-RI-phonon (□), electron- 
electron Coulomb (+) or intervalley electron-phonon (■) scat- 
tering as function of temperature T. It is noted that when 
T > 150 K, only the total SRT and the calculation without 
the intervalley electron-phonon scattering are shown since the 
remaining ones almost coincide with the total SRT. In the cal- 
culation, the electron density N e 



To understand the behavior of the spin relaxation, we 
show in Fig. [T] the temperature dependence of the SRT 
calculated with the electron-electron Coulomb, electron- 
RI-phonon, electron-AC-phonon, long-range electron- 
impurity, or intervalley electron-phonon scattering re- 
moved, separately. We find that at low temperature, 
the intervalley electron-phonon scattering is negligible. 
Therefore, spins relax independently in the two val- 
leys and the SRT is solely determined by the intraval- 
ley spin relaxation channel. Specifically, the long-range 
electron-impurity scattering is dominant, which leads 
to a weak temperature dependence of the SRT. 36 i 50 
With the increase of the temperature, the contribution 
of the long-range electron-impurity scattering becomes 
marginal whereas the intervalley electron-phonon scat- 
tering becomes important. This intervalley scattering 
opens another spin relaxation channel together with the 
Zeeman-like term with opposite effective magnetic fields 
in the two valleys, which is similar to the case in rippled 
single-layer graphene. 50 This Zeeman-like term is given 
by nPi(k) [Eq. (J3|)] approximately. Then the spin relax- 
ation due to the above intervalley spin relaxation channel 
can be approximated by the rate equations^ 



s Mfc (t) - s_ 



fj,k 



(*) 



T v (k) 



= 0, (13) 



where S^fe, u)^ = [^Pi{k)]i and r v (k) represent the in- 
plane spin vectors in each valley, spin precession vec- 
tor along the z-axis and the intervally electron-phonon 



,(*) = 



r v (k) weak scattering 

limit (\^(k)\T v (k) > 1) 

|ft(fc)f a T,(fc) strong scattering 
limit (\^(k)\T v (k) « 1) 



(14) 



Here, r Sy stands for the SRT along the y-axis. Since 
the system we investigate is always in the degenerate 
regime, we take r Sx (k) w t Sj . (kp) approximately with 
kF representing the Fermi wave vector. In our calcula- 
tion, the SRT along the x-sxis is solely determined by the 
intervalley electron-phonon scattering at T ~ 100-400 K. 
At T = 100 K, \p x (k F )\T v {k F ) 160, i.e., the inter- 
valley electron-phonon scattering is in the weak scatter- 
ing limit. As a result, the SRT t S[c = r v (k F ) according 
to Eq. (TT4"]) , which decreases with the enhancement of 
the intervalley scattering as the temperature increases. 
Nevertheless, at T = 400 K, \Pi(k F )\T v (k F ) « 0.04, i.e., 
the intervalley scattering is in the strong scattering limit. 
Then the SRT r Sx = 2/[\p 1 (k F )\ 2 T v (k F )} from Eq. (fli]|, 
which increases with increasing of T. Therefore, the SRT 
first decreases then increases with increasing temperature 
when the intervalley scattering changes from the weak to 
strong scattering limit. The crossover between the weak 
and strong intervalley scattering limit is determined by 
T v {k F )~ l ~ \Pi(k F )\£^ At this crossover point, 



\Pi(k F )\- 



(15) 



Here r Bx rs 315 ps, very close to the value 330 ps shown 
in the figure. It is noted that the electron-electron 
Coulomb scattering is unimportant in the whole temper- 
ature regime. 



2. Electron- density dependence of spin relaxation 

Then wc turn to study the electron-density depen- 
dence. In Fig. [5J the SRT along the re-axis is plotted 
against the electron density at T = 100 (300) K. One ob- 
serves a peak in the electron-density dependence of the 
SRT at T = 100 K whereas the SRT increases with the 
increase of the electron density when T = 300 K. We 
first focus on the case at T = 100 K and show the SRT 
with the intervalley electron-phonon scattering excluded 
in Fig. [5] It is seen that when the electron density is low, 
the intervalley electron-phonon scattering is marginal. 
The SRT is then determined by the intravalley spin relax- 
ation channel and hence increases with the suppression 
of the inhomogeneous broadening^ as the electron den- 
sity increases^ The suppression of the inhomogeneous 
broadening can be understood from the decrease of the 
spin splitting with increasing momentum at high momen- 
tum [see Eqs. HKS])], as shown in Fig. [HI in Appendix |A"1 
This behavior is very different from the case in both semi- 
conductors and single-layer graphenei 36 i 70 When the elec- 
tron density further increases, the intervalley electron- 
phonon scattering becomes important. This intervalley 
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scattering is always in the weak scattering limit with 
\^i(kF)\T v (kp) > 40 and therefore t Sj . = r v (k F ) accord- 
ing to Eq. (fTl]) . As a result, the SRT decreases with the 
increase of the intervalley scattering and hence the in- 
crease of the electron density^ Therefore, in the whole 
electron-density regime, a peak is observed. It is noted 
that this peak is very different from the one predicted 
by Jiang and Wu 7 -^ in semiconductors where the peak 
is attributed to the crossover from the nondegenerate- 
to-degenerate limit. As for the case at T = 300 K, 
the SRT is always dominated by the intervalley electron- 
phonon scattering as shown in Fig. [5] by comparing the 
SRT calculated with and without the intervalley electron- 
phonon scattering. Since this intervalley scattering is in 
the strong scattering limit with \Pi(kp)\T v (kF) < 0.75, 
the SRT t Sx = 2/[\p l {k F )\ 2 T v {k F )] from Eq. This 
SRT increases with increasing electron density since both 
|/3i(fc^)| and the intervalley scattering time r v (kp) de- 
crease with the increase of the electron density 50 ' 55 
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FIG. 2: (Color online) SRT along the a;-axis r Sx as function 
of the electron density N e with x (o) and without □ (a) the 
intervalley electron-phonon scattering at T — 100 (300) K. 



3. Influence of initial spin polarization 

The initial spin polarization dependences of the SRT 
along the £-axis at T = 100 and 300 K are shown 
in Fig. |31 It is seen that the SRT decreases rapidly 
with the increase of initial spin polarization P at T = 
100 K whereas the SRT shows a slight increase with in- 
creasing P at T = 300 K. The rapid decrease of the 
SRT at T = 100 K is very different from the previous 
studies in both semiconductor o 72 ' 76 ' 80 ' 81 and single-layer 
graphene^ where the SRT increases significantly with in- 
creasing initial spin polarization. In the previous studies, 
the increase of the SRT originates from the contribution 
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FIG. 3: (Color online) Total SRT along the z-axis t Sx A (▼) 
and the calculation without the Coulomb HF term o (x) or 
the intervalley electron-phonon scattering • (□) as function 
of the initial spin polarization P at T — 100 (300) K. In the 
calculation, the electron density iV e = 3 x 10 12 cm -2 . 



of the Coulomb HF term.22- This term serves as an effec- 
tive magnetic field along the direction of the spin polar- 
ization, which blocks the spin precession induced by the 
Dresselhaus^ 8 - and/or Rashba^ SOCs and slows down 
the spin relaxation. However, the contribution of the 
Coulomb HF term is marginal to the spin relaxation in 
our investigation, which is shown in Fig. [3] by comparing 
the the SRT calculated with and without the Coulomb 
HF term. Here, the SRT is determined by the inter- 
valley electron-phonon scattering, as shown in Fig. [3] by 
comparing the calculation with and without the interval- 
ley electron-phonon scattering. The intervalley electron- 
phonon scattering transfers electrons between the two 
valleys. The transferred electrons experience opposite ef- 
fective magnetic fields induced by the Zeeman-like term 
in the two valleys, which are not affected by the effective 
magnetic field induced by the Coulomb HF term. More- 
over, the intervalley electron-phonon scattering is in the 
weak scattering limit at low temperature and hence the 
SRT r Sx = r v (k F ) [see Eq. ([I3|l]. The intervalley scatter- 
ing time T v {kp) decreases significantly with the increase 
of the initial spin polarization as the density of states in- 
creases with the polarization. Therefore, the SRT shows 
a rapid decrease with the initial spin polarization at 
T = 100 K. As for the case of T = 300 K, the SRT is also 
determined by the intervalley electron-phonon scattering 
in comparison with the case of T = 100 K. However, in 
contrast to the case of T — 100 K, the intervalley scat- 
tering is in the strong scattering limit. This makes the 
SRT present an opposite trend to the case of T = 100 K. 
A mild increase of the SRT results from the insensitivity 
of T v (kp) to the initial spin polarization at T = 300 K. 
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FIG. 4: (Color online) SRT along the a; (z)-axis r Sa . (r 3 _) as 
function of temperature T. In the calculation, the electron 
density N e = 3 x 10 12 cm" 2 . 



4- Amsotropy of spin relaxation 

We also address the anisotropy of the spin relaxation 
with respect to the spin polarization direction. Due to 
the existence of the out-of-plane effective magnetic field 
[see Eq. |(3J], the in-plane SRTs are identical^ For com- 
parison, we plot the temperature dependence of the SRT 
along the x (z)-axis as the red curve with crosses (blue 
curve with open squares) in Fig. [4] It is seen that in 
comparison with the SRT along the z-axis, the one along 
the x-axis is comparable at low temperature and orders 
of magnitude smaller at high temperature. As men- 
tioned previously, the SRT along the x-axis is strongly 
suppressed by the intervalley spin relaxation channel in- 
duced by the intervalley electron-phonon scattering and 
the Zeeman-like term in the two valleys at high tem- 
perature. However, the contribution of this spin relax- 
ation channel to the SRT along the z-axis is marginal [see 
Eqs. (JTJ2J)] . Therefore, a strong anisotropy between the 
out-of- and in-plane spin relaxations is observed at high 
temperature. 



B. Comparison with experiments 

As mentioned in the introduction, experiments on 
spin relaxation have been carried out in BLG on Si02 
substrat o 37 ! 38 ' 66 or in freely suspended BLG22, very re- 
cently by different groups. The in-plane SRTs in these 
experiments are of the order of 0.01-1 ns, which are 
comparable to our theoretical results at high tempera- 
ture (several hundred picoseconds). This motivates us 
to carry out a detailed comparison with the existing ex- 
periments. In these experiments, the gate voltage, tem- 
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FIG. 5: (Color online) In-plane SRT r 3 as function of tem- 
perature T. + (A) corresponds to the case of the gate 
voltage Vcnp + 20 (60) V in the experiment of Han and 
Kawakami^ x (□) stands for the case of electron density 
N e — 0.75 (1.5) x 10 12 cm -2 in the experiment of Avsar et 
al.& T (■) represents the case of N e = 0.7 (2.2) x 10 12 cm" 2 
in the experiment of Yang et al.^ A denotes the case of 
N e = 10 11 cm" 2 in the experiment of Neumann et a/..— (a) 
Red solid (orange chain) curve stands for our calculation cor- 
responding to the case of Vcnp + 20 (60) V in the absence 
of the short-range scatterers. In the inset, we show the tem- 
perature dependence of the long-range impurity density Ni 
corresponding to the case of Vcnp + 20 V. (b) Grey dashed 
and skyblue double-dotted chain curves are obtained from 
our calculations corresponding to the case of Vcnp + 20 V by 
including the short-range scatterers, with a fixed short-range 
impurity density = 2 x 10 10 cm" 2 and a tunable (shown 
in the inset), respectively. We also show the dependence of t b 
on rid (with the scale on the top of the frame) at T = 50 K 
with other conditions the same as the case of Vcnp + 20 V. 
The short-range potential strength Vq = 6 X 10~ 17 meVm 2 . 
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perature and electron density^ - (tuned by the gate volt- 
age) are given explicitly. However, some parameters such 
as the out-of-plane electric field and the impurity den- 
sity are not available in the experiments. In our calcu- 
lation, the out-of-plane electric field is taken to be the 
ratio of the gate voltage to the thickness of the substrate 
approximatel y 16 i 46 The impurity density is obtained by 
fitting to the mobility or the spin diffusion coefficient 
given in the experiments^ In addition, we assume the 
initial spin polarization along the re-axis without loss of 
generality since the in-plane spin relaxations are identical 
as pointed out previously. 

We first compare the temperature dependence of the 
SRTs obtained from the experiments and calculation 
from our model. In these experiments, the mobilities 
are reported of the order of 10 2 -10 4 cm 2 /(Vs). The 
SRTs from different groups are comparable as shown 
in Fig. However, the SRT from the experiment of 
Han and Kawakami^l decreases with increasing temper- 
ature whereas the ones from the experiments of Yang 
et alf& and Avsar et al^ show a marginal tempera- 
ture dependence. In addition, Han and Kawakami^ 7 - re- 
ported the longest SRT among these experiments, in- 
dicating that their experimental data may contain less 
extrinsic effects. Motivated by this, we show our com- 
parison with the experiment of Han and Kawakami 3 - 7 - in 
Fig. EJa). The experimental data with the gate voltage 
Vcnp + 20 (60) V corresponding to the electron density 
N e = 1.47 (4.4) x 10 12 cm~ 2 , is labeled as + (A). By 
fitting to the spin diffusion coefficient given in the ex- 
periment, we obtain the long-range impurity density Ni 
shown in the inset as the curve with •. In contrast to 
the long-range impurity scattering, the contribution of 
other scatterings to the spin diffusion coefficient (mobil- 
ity) is marginal. Then we calculate the SRT correspond- 
ing to Vcnp + 20 (60) V shown as the red solid (orange 
chain) curve. A crossover is observed in the tempera- 
ture dependence of the SRT, similar to the one shown 
in Sec. IIIA. It is seen that the SRT from our calcula- 
tion is comparable to the experimental data at high tem- 
perature. This suggests that the intervalley spin relax- 
ation channel induced by the intervalley electron-phonon 
scattering and the Zeeman-like term in the two valleys 
addressed in Sec. IIIA, plays an important role at high 
temperature. Nevertheless, the SRTs are still larger than 
the experimental data with an observable difference at 
high temperature. Additionally, at low temperature, our 
results are orders of magnitude larger than the experi- 
mental ones since the intervalley electron-phonon scat- 
tering is marginal at low temperature and the SRT is 
determined by the long-range electron-impurity scatter- 
ing in our calculation. Due to the existence of all these 
discrepancies, one has to take into account extrinsic ef- 
fects such as adatomS ) 24 ' 42 " — cur vat ur e 1 1 6 ' 48 ~ ^£ substrate 
effects^-^ and also contacts . 30 ' 39 ' 40 ' 53 In the above 
experiments on BLG, the spin polarized electrons are 
all injected from ferromagnetic contacts! 37 ' 38 ' 66 ' 68 As re- 
ported, the ferromagnetic electrodes generate stray fields, 



which cause spin precession around the contact) 5 - 3 - As a 
result, the stray fields may have a strong influence on 
the observed experimental results. Additionally, electron 
scattering by the tunneling contacts may also be respon- 
sible for the available experimental datai 30 ' 39 ' 40 

In addition to the above extrinsic effects, the short- 
range scatterers may also be an extrinsic sourcei 9 ' 92 ' 93 
Here, we discuss the contribution of short-range 
scatterers 9 - in detail. The scattering term induced by the 
short-range scatterers is given by 

<9 t /3pk|sR = ~Tm d Vu ^2 l^k Vk'| 2( 5(Vk' - e ^k) 

x (ft lk -/Vk') + H - c -> (16) 

with rid and Vq denoting the short-range impurity den- 
sity and the constant short-range potential strength, 
respectively— This scattering term contributes to both 
the intra- and inter- valley scatterings. The intervalley 
scattering opens an intervalley spin relaxation channel 
together with the Zeeman-like term in the two valleys [see 
Eq. ((3J)], similar to the intervalley electron-phonon scat- 
tering. However, in contrast to the intervalley electron- 
phonon scattering, the intervalley scattering induced by 
the short-range scatterers is insensitive to the tempera- 
ture and may play an important role in the in-plane spin 
relaxation especially at low temperature. Here, we choose 
a fixed short-range impurity density nj = 2x 10 10 cm~ 2 
and potential strength Vb = 6x 10~ 17 meVm 2 . It is noted 
that although with this additional short-range scatter- 
ing included, the mobility is still dominated by the long- 
range electron-impurity scattering since the contribution 
of the long-range electron-impurity scattering is about 
four orders of magnitude larger than that of the short- 
range one. We calculate the temperature dependence 
of the SRT corresponding to the case of Vcnp + 20 V 
in the experiment of Han and Kawakami 37 by including 
the additional short-range scattering. The result is plot- 
ted as the grey dashed curve in Fig. [Sfb). A crossover 
is also observed at high temperature due to the contri- 
bution of the intervalley electron-phonon scattering as 
pointed out previously. By comparing our calculation 
with the experimental data labeled as +, we find that our 
result becomes comparable at low temperature. Further- 
more, by fitting to the experimental data with tunable 
short-range impurity density shown as the curve with o in 
the inset, our result labeled as the skyblue double-dotted 
chain curve agrees fairly well with the experimental data 
in the temperature regime lower than the crossover. As 
for the temperature higher than the crossover, our result 
becomes even larger than the experimental one once the 
short-range scattering is included. This is because that 
the intervalley electron-phonon scattering is in the strong 
scattering limit. With the inclusion of the short-range 
scattering, the intervalley scattering becomes stronger 
and the SRT becomes larger. This discrepancy may be 
attributed to other extrinsic effects. 

As mentioned above, the SRT can be tuned by the 
short-range impurity density. We show the short-range 
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impurity density dependence of the in-plane SRT at 
T = 50 K as the purple dotted curve in Fig. [5jb), with 
all other conditions remaining the same as the case of 
Vcnp + 20 V. A valley is also observed, similar to the 
case of the temperature dependence discussed above. At 
low temperature, the intervalley electron-phonon scat- 
tering is negligible. The SRT is then determined by the 
intervalley short-range scattering. The valley originates 
from the crossover between the weak and strong inter- 
valley short-range scattering limit by tuning the short- 
range impurity density. At the crossover point, the SRT 
t s s» 232 ps according to Eq. (fT5|) . close to the value 
252 ps shown in the figure. This crossover point gives the 
minimum of the SRT and can be taken as a criterion for 
determining the possible contribution of the short-range 
scattering. Specifically, the crossover point r s w 316, 
258, 232, 283, and 249 ps corresponding to the case of 
the experimental data shown as ▲, x, □, T and I, re- 
spectively We find that both the experimental data of 
Han and Kawakami 37 and Avsar et al£2- are larger than 
the crossover point whereas the one of Yang et al£^ is 
smaller. This suggests the experimental data of Han and 
Kawakami 3 ^ and Avsar et alM- can be explained by the 
short-range scatterers. As for the experiment of Yang et 
al.r^- some other extrinsic effects have to be considered. 




Yang et al. 
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FIG. 6: (Color online) In-plane SRT t s as function of electron 
density. The dotted curves with symbols correspond to the 
experimental data of Han and Kawakami (HK)' 7 at T = 20 K 
and Avsar et a/.— and Yang et aZ.— at T = 5 K. The solid 
curves with symbols correspond to the experimental data at 
room temperature of Han and Kawakami,— Avsar et a/.— and 
Yang et aZ..— Red solid (blue chain) curve is from our calcu- 
lation corresponding to the case of room (low) temperature in 
the experiment of Han and Kawakami.— It is noted that the 
blue chain curve is calculated with the short-range scatterers 
with the short-range impurity density rid = 10 10 cm -2 and 
potential strength Vo = 6 x 10~ 17 meVm 2 . 

Then we turn to compare the electron density depen- 
dence of the SRTs from the experiments and our calcu- 



lation. Among the experiments, Han and Kawakami 37 
measured the longest SRT at low temperature (up to 
several nanoseconds) shown as the dotted curve with A 
in Fig. [fpj This SRT is about one order of magnitude 
larger than the ones from other experimental cases. In 
addition, the SRTs from the experiments show different 
dependences on electron density. Specifically, at low tem- 
perature, both Avsar et al^ and Yang et al£H show a 
decrease of the SRT with the increase of the electron 
density whereas a mild peak is observed by Han and 
Kawakami i^I At room temperature, the SRTs from Avsar 
et al£2- and Yang et al^ increase with increasing elec- 
tron density whereas Han and Kawakami 37 reported an 
insensitive density dependence of the SRT. 

We first show our comparison with the experimen- 
tal data at low temperature. As shown in Fig. [5ja), 
our result is orders of magnitude larger than the ex- 
perimental data of Han and Kawakami^ 7 - at low tem- 
perature without short-range impurities. However, with 
the inclusion of the short-range scatterers, our result 
shows good agreement with the experimental data. Here, 
with the short-range scatterers explicitly included, we 
also show our comparison with the experimental data 
of Han and Kawakami^ 7 . With a fixed short-range im- 
purity density — 10 10 cm' 2 and potential strength 
Vo = 6 x 10~ 17 meVm 2 , the electron density depen- 
dence of the SRT from our calculation is shown as blue 
chain curve in Fig. [SJ We find that our result marginally 
agrees with the experimental data shown as the dotted 
curve with ▲. Better agreement can be obtained with 
tunable short-range impurity density. In addition, we 
also find that the experimental data of Avsar et al^ can 
be explained by the short-range scatterers. However, to 
account for the experimental data of Yang et al.^ one 
has to include other extrinsic effects. This is because the 
SRT of Yang et al£^ is smaller than the crossover point in 
the short-range impurity density dependence of the SRT 
(230-320 ps), which is taken as a criterion for the possible 
contribution of the short-range scattering as mentioned 
previously. 

We then address our comparison with the experimen- 
tal data at room temperature. As shown in Fig. [5ja), 
our result is larger than the experimental data of Han 
and Kawakami^ 7 . with an observable difference at room 
temperature without short-range impurities. With the 
short-range scatterers included, the SRT from our calcu- 
lation becomes even larger than the experimental data 
as mentioned previously, indicating that other extrinsic 
effects have to be taken into account. Here, we show our 
comparison with the experiment of Han and Kawakami 37 
without extrinsic effects. The experimental data and our 
result are shown as the solid curve with • and the red 
solid curve in Fig. [JoJ respectively. We find that the SRT 
in our calculation is about three times as large as the 
experimental one in low electron density regime and one 
order of magnitude larger than the experimental one in 
high density regime. 
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IV. SUMMARY 

In summary, we have investigated the electron spin re- 
laxation due to the DP mechanism in BLG with only 
the lowest conduction band being relevant. The SOC 
of the lowest conduction band is constructed from the 
symmetry group analysis with the parameters obtained 
by fitting to the numerical calculation according to the 
latest report by Konschuh et al£^ from first principles 
with both the intrinsic and extrinsic SOC terms in the 
pseudospin space included. We find that the magni- 
tudes of both the out-of- and in-plane components of 
the SOC decrease with increasing momentum at large 
momentum, indicating a suppression of the inhomoge- 
neous broadening with the increase of the momentum. 
This is different from the case in both semiconductors 
and single-layer graphene. Additionally, the leading term 
of the out-of-plane component of the SOC serves as a 
Zeeman-like term with opposite effective magnetic fields 
in the two valleys, which is similar to the case in rip- 
pled single-layer graphene. This Zeeman-like term, to- 
gether with the intervalley electron-phonon and/or pos- 
sible intervalley short-range scatterings, opens an inter- 
valley spin relaxation channel, which has not been re- 
ported in the literature in BLG. The intervalley electron- 
phonon scattering is derived by using the tight-binding 
model. In addition to the intervalley electron-phonon 
and short-range scatterings, we also include the long- 
range electron-impurity, electron-electron Coulomb and 
intravalley electron-phonon scatterings to calculate the 
SRT. We find that the in-plane SRT is strongly sup- 
pressed by the intervalley electron-phonon scattering at 
high temperature. In contrast to the intervalley electron- 
phonon scattering, the intervalley short-range scattering 
is insensitive to the temperature and plays an important 
role in the in-plane spin relaxation especially at low tem- 
perature. 

In the absence of short-range scatterers, a marked non- 
monotonic temperature dependence of the in-plane SRT 
is predicted with a minimum SRT down to several hun- 
dred picoseconds. This nonmonotonic behavior origi- 
nates from the crossover between the weak and strong in- 
tervalley electron-phonon scattering. Moreover, we pre- 
dict a peak in the electron density dependence of the 
in-plane SRT at low temperature, which is very different 
from the one in semiconductors. At high temperature, 
the in-plane SRT increases monotonically with increas- 
ing density. We also find that the in-plane SRT decreases 
rapidly with the increase of the initial spin polarization at 
low temperature. This is very different from the previous 
studies in both semiconductors and single-layer graphene 
where the SRT increases significantly with increasing ini- 
tial spin polarization. The physics is understood that 
the spin relaxation time is determined by the interval- 
ley electron-phonon scattering, which transfers electrons 
between the two valleys. The transferred electrons expe- 
rience opposite effective magnetic fields in the two val- 
leys, which are not affected by the Coulomb HF term. 



In addition, a strong anisotropy between the out-of- and 
in-plane spin relaxations is also addressed at high temper- 
ature where the out-of-plane SRT is about two orders of 
magnitude larger than the in-plane one. As for low tem- 
perature, the out-of- and in-plane spin relaxation times 
are comparable. 

We also show our comparison with the existing exper- 
iments of Han and Kawakami, 37 Avsar et al£& and Yang 
et al.£^ We find that without intervalley scattering, the 
SRT is orders of magnitude larger than the experimental 
data in the whole temperature regime. With the inter- 
valley electron-phonon scattering explicitly included, the 
SRT from our calculation becomes comparable to the ex- 
perimental data at high temperature but still orders of 
magnitude larger than the experimental data at low tem- 
perature. In addition, a crossover in the temperature de- 
pendence of the in-plane SRT is shown, which also results 
from the crossover between the weak and strong interval- 
ley electron-phonon scattering. With the inclusion of the 
short-range scatterers, our result agrees fairly well with 
the experimental data in the temperature regime lower 
than the crossover. As for the temperature higher than 
the crossover, other extrinsic effects have to be included. 
Moreover, a crossover point is also shown in the short- 
range impurity density dependence of the in-plane SRT, 
which gives a minimum for determining the possible con- 
tribution of the short-range scattering. The experimental 
SRT larger than the minimum can be explained by the 
possible short-range scatterers whereas other extrinsic ef- 
fects have to be considered for experimental data smaller 
than the minimum. 
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Appendix A: SOC OF THE LOWEST 
CONDUCTION BAND 

The SOC term of the lowest conduction band near the 
Dirac points can be described by a 2 x 2 matrix 

ttSO _ f n /ill n pX2 \ / 4 i\ 
n efin — rrSO * ttSO I 
V n fj,12 n fj,22 J 

in the basis i/)£ s = Vk ® \s) where s = and ip£ = 

^ iC f(k)*f(k) (i = Ai, Bi, A 2 , B 2 ). The coefficients 
cf (k) can be obtained by exactly diagonalizing the 4x4 
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FIG. 7: (Color online) (a) a\, (b) a\ and (c) a? 3 as function 
of the out-of-plane electric field E z with j — 1-5. It is noted 
that a}, of, af, af and af (i = 1-3) are dimensionless. 



effective Hamiltonian given in Eq. ([5|). Here, 



ttSO 

"■till 



rrSO _ 



ttSO 



X;<*(k)^(k)<^ t |Hso|^, t ), (A2) 
^ C r(k)^(k)<^|ff S o|^}, (A3) 
£ C f(k)c?(k)^ t |ff S o|^>, (A4) 



with (y^Wsol^J 0i,2 = {t4}) being the spin- 
orbit matrix elements near the Dirac points in the pseu- 
dospin space given by Konschuh et aL— Specifically, 



H 



so 
tin 

so 

p22 



= ^\ n {\4\-\^\ 2 )+^i2{\c>t\ z -\4\ 2 ) (A5) 

(A6) 



rrSO 
--"/ill 



frSO 



[ C 1* C 3 (^4) + C2*Ci ( — iAo) + C^*^ — £A 4 ) 



+ c£ c£(— iAo) + c%(-i\' 4 ) + c% C3 (-iA 3 ) 



J 1 * J 1 



c^(iA^)] 



(A7) 



An and A12 are the strengths of single-layer-like intrinsic 
SOC; Ao and Aq represent single-layer-like extrinsic spin- 
orbit strengths; A3, A4 and A 4 stand for interlayer spin- 
orbit parameters^ Then the effective magnetic field of 
the SOC ft^(k) in Eq. @ is given by fi£(k) = 2ReH^ 2 , 
fi£(k) = -2I mj ffSO and QM(k) = 2H^ V It is noted that 
the contribution of other energy bands to the SOC of the 
lowest conduction band is marginal in our calculation. 

As reported by Konschuh et al.^ the small group of 
the Dirac points is C3 in the presence of an out-of-plane 
electric field. Following from the symmetry group C^r^ 
we derive an analytical form of the SOC near the Dirac 
points [see Eqs. (p][3])]. In these equations, the coeffi- 
cients ai(k) (i — 1-3) and (i — 1-2) are obtained by 
fitting to the numerical calculation from Eqs. (|A5|) and 
(|A7|) using the Pade approximation^ Specifically, the 
coefficients a^(fc) (i = 1-3) in Eqs. ([TJ2]) read 



ctfi(fc) = Anafc 



a 7 



afak 



1 + afak 



fa 2 k 2 



(A8) 



The coefficients /3,(fc) (i = 1-2) in Eq. © are given by 

I3\ + Rah 



fh{k) = An a/c 



1 + (3fak + (3fa 2 k 2 ' 
ft + /3 2 ak 



1 + Plak + fila 2 k 2 ' 



(A9) 
(A10) 



1,3 



Here, the coefficient 0® = —24 fieV, which is indepen- 
dent of the applied electric field. The electric field de- 
pendences of a\ (i — 1-3, j — 1-5) and (3f (i = 1-2, 
j = 1-4) are plotted in Figs. [7] and [5J respectively. It is 
noted that the Pade approximation^ can give a precise 
description of (i = 1-3) and (3i(k) (i = 1-2) at large 

momenta near the Dirac points, which can be populated 
and play an important role in spin relaxation according 
to the experimental conditions (heavily doped) i 37 i 38 ' 66 i 68 
With only a simple linear fc-order approximation, aj(fc) 
and Pi{k) agree with the numerical results only at very 
small momentum near the Dirac points. 

We also show in Fig. [9] the spin splitting of the lowest 
conduction band calculated with analytical form of the 
SOC (blue dashed curve) and with the explicit numer- 
ical one (red solid curve). We find that the analytical 
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FIG. 8: (Color online) (a) PI and (b) P 3 2 as function of the 
out-of-plane electric field E z with j — 1-4. It is noted that 
Pi, /3f, fif and Pi (i = 1,2) are dimensionless. 



result agrees fairly well with the numerical one at large 
momentum. In addition, the spin splitting decreases 
with increasing momentum at large momentum, which 
is very different from the case in both semiconductors 70 
and single-layer graphene^ 



Appendix B: SCATTERING MATRIX ELEMENTS 



The electron-electron Coulomb scattering matrix ele- 
ment is given by |V^ k _ q | 2 . The screened Coulomb po- 
tential V£ k = V^/eC^e/ik - e M k-q) with e(q, e Mk - 



) = 1 - Vq (0) Il(q, e Mk - e M k- q ) being the screen- 



ing under the random phase approximation^ V^ ' = 
2Trvpr s /q is the two-dimensional bare Coulomb poten- 
tial with vf and r s being the Fermi velocity in single- 
layer graphene^ 7 . and the dimensionless Wigner-Seitz 



FIG. 9: (Color online) Spin splitting of the lowest conduc- 
tion band AE near the K point as function of the momen- 
tum k x under a typical electric field E z — 0.14 V/nm. Solid 
(dashed) curve represents the numerical (analytical) result. 
|K] =47r/(3a). 



radius^ — respectively. II(q, ui) is given by 



36.99-101 



n(q,w) 



E\rpliw' 
l J kk+q 

^ikvv' s 



h 



ks 



J k+q.s 



+ uj + i0+ 



(Bl) 



f/)^" Vk+q- The long-range electron- 



where T^ q 

impurity scattering matrix element |^kk-ql 2 = 

Z i\ V kk~q\ 2e ~ 2qd ' 111 wllich Z i and d stand f° r thc im- 
purity charge number and effective distance of the im- 
purity layer to the BLG sheet 44,96-98,102 It ig noted 

that when we calculate the electron-electron Coulomb 
and long-range electron-impurity scatterings, we take 
the distance between the graphene layers to be zero 
approximately. 10 = It is further noted that the Coulomb 
potential in the HF term is in the static screening limit, 

i.e,^ q = ^Mq,0). 

Then we turn to the electron-phonon scattering ma- 
trix elements. For the intravalley electron- AC-phonon 
scattering, the scattering matrix element \M^ k , k , \ 2 = 

^AC^/( 2 P w piO / kk'<W where q = |k - k'|; D AC and 
v p h represent the deformation potential— and acoustic 
phonon velocity ; 104 ' 105 respectively. 

For the intravalley electron-RI-phonon scattering, 

KkVk'l 2 = S^r*££/kV with q s = Ar s k F be- 
ing the Thomas-Fermi screening length^ For Si02 sub- 
strate, the energy spectra of the two remote phonon 
modes are denoted as ujf~ l and w^ 1 , respectively; g\ 
and (72 are the corresponding dimensionless coupling 
parameters 

The matrix elements for the intravalley electron-OP- 
phonon scattering (\M° k , k , | 2 ) are derived using the 
tight-binding model according to the arXiv version of 
Ref. |63. The intravalley electron-OP-phonon scattering 
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includes both the electron in-plane (\Mj£ , k , | 2 ) and out- 

of-plane (\Mffc„, k t\ 2 ) OP phonon scatterings. The elec- 
tron in-plane OP phonon scattering matrix element is 
given by 



with ilzo being the energy of the out-of-plane OP phonon 
mode.— 



9h 2 



{lV'k t [7o(sin6 l qCT r j 3 + cos^qcr^ 3 ) Appendix C: RELEVANT DIRAC MATRICES 



— aj' 3 (sin 6> q 7Q7Q — i sin q 7o + i cos # q 7o + cos 6> q 
x lhlh)/(2V3k)m 2 + \^H(^ Vd - cos0 q 
x + af'Scose^l - sin ^71,7^/(^4)]^ I 2 
+ l^k f [7o(cos Vd 3 " smO^ 3 ) - a l3 (cos q7 fj7D 

- tcosB^l -isinO^h - sin ^7^/(273/3)]^ | : 
+ iV^N^cosflqCT^ 1 + sin<V° 2 ) - aitfl (i sin6» q 



-I- cos 



(B2) 



where Olt is the energy of the in-plane OP phonon 
modes j 83 ' 84 ^3 is the bond length corresponding to the 
interlayer hopping 73; and o^ 2 and c^ 3 are 4x4 Dirac 
matrices given in Appendix [Cj 86 Here, # q is the polar 
angle of the momentum q. The electron out-of-plane OP 
phonon scattering matrix element reads 



lAf, 



zo 

/xk,/j'k' I 



Zpilzo 



The relevant Dirac matrices^ used in the electron- 
phonon scattering are 



7d = 



,01 



^02 

D 



I 

1 



a x 

-<7, 



a 







7d = 



,13 



^23 
,cr D 



-/ 
I 



(Jy 

cr„ 



(T^ 

CTj; 



(B3) with / being 2x2 unit matrix. 



(CI) 

,(C2) 
(C3) 

(C4) 
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